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Chemical  reactions  between  ions  and  electrons  of  a  translating,  high  velocity  field  re¬ 
versed  plasmoid,  including  charge  exchange  and  electron  impact  ionization,  are  examined 
in  an  adiabatic  heat  bath  relaxation.  Importance  of  different  reaction  mechanisms  is  ana¬ 
lyzed  for  conditions  typical  for  FRC  thrusters,  with  gas  and  plasma  densities  on  the  order 
of  1018  molecule/m3,  plasma  temperatures  between  5  and  50  eV,  and  relative  velocities 
up  to  30  km/s.  Then,  an  implicit  Particle-In-Cell  code,  Celeste3D,  extended  to  include 
plasma-neutral  and  neutral-neutral  elastic  and  inelastic  processes  and  Coulomb  collisions, 
is  used  to  model  a  planar  flow  between  an  FRC  plasmoid  initially  at  a  Schmid-Burgk  equi¬ 
librium,  and  neutral  gas  with  varying  relative  velocities,  plasma  temperatures,  and  neutral 
gas  densities.  The  contribution  of  charge  exchange  and  electron  impact  ionization  reactions 
is  analyzed,  and  the  sensitivity  of  the  entrainment  efficiency  to  the  plasmoid  translation 
velocity  and  neutral  density  is  evaluated. 


I.  Introduction 

One  of  the  most  promising,  and  at  the  same  time  challenging,  propulsion  technologies  is  high-power 
electric  propulsion.  The  first  challenge  is  associated  with  the  necessity  to  operate  a  thruster  at  very  high 
efficiency,  otherwise  the  system  limitations  due  to  heat  rejection  become  insurmountable.  It  is  clear  that 
high  power  can  be  efficiently  delivered  into  a  plasma  by  raising  its  temperature,  i.e.  its  specific  energy; 
consequently,  the  specific  impulse  can  take  large  values,  even  in  the  excess  of  10,000  sec.  However,  for  a  fixed 
power  available,  this  also  translates  into  small  thrust  and  acceleration.  While  such  a  high  Isp  may  be  highly 
desirable  for  deep  space  missions,  the  low-Earth  orbit /Geosynchronous  Earth  orbit  environments  in  which 
the  US  Air  Force  operates  is  mostly  characterized  by  two  competing  requirements,  i.e.  (1)  minimization  of 
propellant  and  system  mass  (high  Isp)  and  (2)  minimization  of  the  time  required  for  achieving  the  orbital 
maneuver  (high  thrust).  Following  the  “goldilocks”  principle,  one  must  find  the  optimal  range  of  Isp  which, 
depending  on  the  mission,  is  usually  between  2-5  ksec.  The  second  challenge  is  therefore  to  design  an  electric 
propulsion  (EP)  system  which  can  operate  in  that  regime  at  high  power  (100-200  kWe).  Finally,  it  would 
be  beneficial  to  to  operate  in  a  dual  mode,  i.e.  be  able  to  tune  the  regime  between  higher  thrust/lower  Isp 
and  lower  thrust/higher  Isp,  at  constant  power  and  efficiency,  thus  presenting  a  third  challenge. 

Plasma  at  high  temperature  inevitably  becomes  highly  ionized;  this  is  very  beneficial  for  coupling  with 
external  electro-magnetic  fields  and  high-efficiency  thrust  generation.  Typical  losses  associated  with  the 
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plasma  in  these  conditions  include:  (a)  radiation;  (b)  effective  ionization  cost;  (c)  loss  of  confinement.  The 
amount  of  energy  lost  to  radiation  depends  strongly  on  the  conditions,  but  generally  speaking,  the  line 
transitions  are  self-absorbed  while  the  plasma  is  closer  to  being  optically  thin  for  continuum  radiation.  The 
integrated  bremstrahlung  emission  is  only  weakly  proportional  to  temperature  and  the  radiative  capture 
rate  (free-bound  transitions)  decreases  at  high  temperature;  nevertheless,  the  rate  of  radiative  cooling  can 
be  problematic  for  high -Z  plasma  (due  to  a  Z 2  dependence)  and  in  radiative  non-equilibrium  conditions 
(volumetric  emission).1  The  effective  ionization  accounts  for  the  amount  of  energy  required  to  produce  an 
ion  by  electron  impact,  accounting  for  all  other  couplings  which  can  reduce  the  electron  energy;  these  include 
elastic  collisions  as  well  as  all  electronic  excitations,  and  vibrational  excitations  in  the  case  of  molecular 
plasma.  From  dependence  with  respect  to  electron  energy  of  the  cross-sections  for  these  various  processes,  the 
most  efficient  regime  is  obtained  for  electron  temperatures  of  the  order  of  50  eV.  Finally,  particle  confinement 
strongly  depends  on  applied  magnetic  fields,  since  the  cross-field  diffusion  scales  approximately2  as  l/B. 

Therefore,  the  optimal  conditions  are  achieved  for  low -Z  plasma  at  high  (50-100  eV)  temperature  and 
in  strong  magnetic  fields.  These  conditions  are  typically  obtained  in  Field  Reversed  Configuration  (FRC) 
plasma.  The  FRC  is  a  self-organized  magnetized  plasma  structure  in  the  shape  of  a  highly  compact  toroid. 
The  magnetic  field  is  mostly  in  the  poloidal  direction,  generated  by  internal  (toroidal)  currents.  The  ratio 
of  plasma  pressure  to  magnetic  pressure  (/3)  is  close  to  unity,  i.e.  the  highest  plasma  density  that  can  be 
attained  for  given  external  magnets.  The  poloidal  field  also  contributes  to  the  particle  confinement.  Starting 
from  a  background  uniform  plasma  at  a  constant  axial  (bias)  field,  the  FRC  can  be  formed  by  pulsing 
external  coils  and  reversing  the  applied  field,  inducing  currents  at  the  plasma  boundary  and  pinching  the 
plasma  at  both  ends  (see  1,  left).  The  initial  bias  field  is  trapped  inside  the  plasma  and  forced  to  reconnect 
at  the  end  points  (separatrix),  creating  an  elongated  toroidal  shape  (see  1,  right). 


Figure  1.  Schematics  of  basic  FRC  formation  by  field  reversal  (left)  and  schematics  of  an  FRC  plasmoid  after 
formation,  where  the  crosses  indicate  the  toroidal  plasma  current  (right). 


In  contrast  to  some  other  high-/sp  plasma  propulsion  concepts,3  the  FRC  is  completely  magnetically 
insulated  from  the  external  field;  in  other  words,  the  plasma  is  not  tied  to  an  external  field  line,  and  the  FRC 
can  readily  detach  from  the  confining  external  field.  It  can  also  be  translated  and  accelerated  by  applying  a 
gradient  of  magnetic  pressure  using  pulsed  external  coils.  The  FRC  can  therefore  be  efficiently  accelerated 
to  provide  thrust,  operates  at  a  temperature  that  is  optimal  for  ionization  and  is  well  confined.  The  basic 
concept  of  operations  has  been  demonstrated4  and  more  recent  research  has  led  to  further  optimization  of  the 
formation  process,5  while  methods  for  increased  efficiency  through  energy  recovery  in  the  electrical  circuits 
are  currently  being  developed  under  an  AFOSR-STTR  program. 

One  of  the  latest  design  iterations6  provides  plasma  velocities  in  the  10-40  km/sec  range,  a  desirable 
regime  for  USAF  applications.  However,  power  and  mass  utilization  efficiencies  still  remain  to  be  determined 
more  precisely  and  thrust  augmentation  is  still  desired.  Since  the  plasma  density  in  the  FRC  is  close  to 
optimal  (/3  ~  1),  this  could  be  achieved  in  two  ways:  (a)  increasing  the  molecular  mass  of  the  propellant;  (b) 
increasing  the  effective  mass  by  gas  entrainment.  The  first  approach  is  by  far  the  easiest  to  achieve  and  one 
would  presume  that  the  common  propellant  of  choice  for  EP  (xenon),  having  the  largest  mass  of  any  (stable) 
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noble  gas  would  be  an  ideal  candidate.  However,  as  previously  mentioned,  a  high-Z  plasma  has  higher  rates 
of  radiative  and  collisional  (effective  ionization  cost)  energy  losses,  thus  reducing  the  efficiency.  The  second 
approach  is  tentatively  more  attractive  but  more  difficult  to  achieve.  In  this  case,  the  FRC  propagates  into 
an  ambient  gas  and  captures  or  pushes  this  neutral  gas  while  being  accelerated  by  the  external  coils.  In 
this  fashion,  the  system  is  able  to  exert  a  force  on  a  composite  object  with  a  higher  mass,  while  the  plasma 
itself  is  kept  at  optimal  physical  conditions  for  high  efficiency.  The  dynamics  of  the  entrainment,  momentum 
coupling  and  mixing  processes  at  the  FRC/gas  interface  are  critical  to  the  efficiency  of  this  concept.  It 
should  be  emphasized  here  that  this  process  is  also  key  to  determining  the  true  mass  utilization  efficiency  of 
the  FRC  thruster  by  itself:  there  is  always  some  level  of  residual,  ambient  neutral  gas  in  the  chamber  as  a 
consequence  of  gas  injection,  whether  from  the  previous  pulse  or  the  current  one.  By  deliberately  introducing 
more  gas  to  be  acted  upon  by  the  FRC  formed  at  the  end  of  the  injection  cycle,  one  can  thus  potentially 
improve  the  thrust  density  while  remaining  at  peak  efficiency. 

In  this  work,  the  process  of  neutral  entrainment  of  an  FRC  plasmoid  is  studied  numerically,  with  the 
main  focus  on  the  analysis  of  the  relative  importance  of  the  electron  impact  ionization  and  charge  exchange 
reactions  between  the  neutral  and  charged  particles  First,  the  balance  between  ionization  and  charge  ex¬ 
change  reaction  rates  is  examined  for  various  gases  and  temperatures.  Then,  a  study  of  heat  bath  relaxation 
is  performed  and  the  impact  of  Coulomb  collisions  is  clarified.  Finally,  a  Celeste3D  particle-in-cell  compu¬ 
tational  tool,'  extended  to  include  the  interactions  between  the  charged  and  neutral  particles  and  neutral 
particle  transport,  as  well  as  Coulomb  interactions  between  charged  particles,  is  used  in  simulations  of  an 
FRC  plasmoid  and  neutral  gas  interaction. 


II.  Rates  of  reactions  involving  charged  particles 


According  to,4,8  the  plasmoid  temperature  is  expected  to  be  less  than  20  eV,  and  in  some  cases  may  cool 
to  as  low  as  5  eV.  Let  us  examine  the  rates  of  the  different  processes  occurring  in  plasma  at  the  temperature 
range  of  interest.  The  following  processes  are  considered: 


1.  Electron  impact  ionization  (Eli,  A +E  — >•  A+  +  e  +  e).  The  electron  impact  ionization  rates  keu 
are  obtained  by  the  integration  of  the  translational  energy-dependent  cross  sections  aeu  tabulated  in 
SIGLO  database9  as 

keii=  J  9<7eii{g)fe{g) dg,  (1) 

where  g  is  the  relative  collision  velocity  and  fe  is  the  Maxwellian  distribution  function. 

2.  Single  charge  exchange  (SCX,  A+  +A  — >  A+A+)  For  charge  exchange  reaction  cross  section,  an 
expression  from10  is  used,  written  in  CCS  units  as 


n  Sir E±JL.\  f> 
u.o7r«o  E_  log  lyuuy  e  mJ  e  ^ 


7r  yj  ae2/2e 

integrated  to  obtain  the  reaction  rate  similar  to  Eq.  1.  Here, 

2ae2  E2 


e 

e  <  e* 


a40  Ry2 


lo§  ^100\/?^- 


e  me 


(2) 


(3) 


is  a  transcendental  equation  for  e*,  ao  is  the  Bohr  radius,  £)  is  the  ionization  energy,  m  and  y  are 
particle  mass  and  reduced  mass,  respectively,  and  a  is  the  polarizability. 

3.  Recombination,  including  the  radiative  recombination  in  three-body  collisions  (A+  +  e  +  e— lA  +  e) 
and  the  radiative  photorecombination  (A+  +  e  +  e  — 1  A  +  e).  The  expressions  from11  are  used,  which 
are,  for  the  photorecombination,  kr i  =  2.7  x  10-19Te-0-75  (nr3/s),  and  for  the  three-body  recombination 
kr 2  =  8.75  x  10_39TfT4'5ne  (m3/s)  For  the  recombination  reactions,  an  ion  and  electron  number  density 
of  1018  nr-3  is  assumed,  similar  to  that  observed  in  experiments  of  Ref.4 


Charge  exchange  reactions  are  beneficial  for  neutral  entrainment,  because  they  produce  a  fast  neutral 
particle  that  contributes  to  thrust,  and  a  slow  charged  particle  that  can  be  quickly  and  efficiently  accelerated 
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by  the  electric  field.  On  the  contrary,  the  electron  impact  ionization,  being  highly  endothermic  reaction, 
takes  away  energy  from  the  system  and  therefore  lowers  efficiency  of  the  thruster. 

The  calculated  reaction  rates  of  these  processes  as  function  of  temperature  for  light  gas,  helium  are 
shown  in  Fig.  2  (left).  It  can  be  concluded  that  at  those  temperatures  and  a  number  density  on  the  order 
of  1018  m— 3,  the  recombination  is  not  important.  The  ionization  reaction  rate  is  lower  than  the  charge 
exchange  rate,  which  may  be  an  indication  that  the  entrainment  process  in  a  helium-based  FRC  thruster 
may  be  expected  to  be  fairly  efficient,  especially  for  lower  plasma  temperatures.  Note  however  that  for 
temperatures  below  10  eV  and  a  typical  plasmoid-neutral  interaction  time  of  50  p, s,  the  fraction  of  ions  that 
participate  in  a  charge  exchange  reactions  may  be  as  low  as  10%  when  the  neutral  density  does  not  exceed 
1018  m— 3.  A  higher  neutral  density  will  therefore  be  required  to  increase  the  number  of  charge  exchange 
reactions  and  thus  improve  the  efficiency  of  the  neutral  entrainment  process.  It  can  be  seen  that  the  ratio  of 
Eli  rate  and  SCX  rate  for  helium  is  favorable  for  neutral  entrainment,  however,  the  thrust  of  a  helium-based 
thruster  is  relatively  low,  and  using  heavier  species  may  be  attractive. 

One  of  potential  candidates  could  be  nitrogen,  especially  keeping  in  mind  its  air-breathing  potential.  The 
charge  exchange,  ionization,  and  recombination  rates  for  nitrogen  are  shown  in  Fig.  2  (right).  Note  that  for 
recombination  rate,  the  dissociative  recombination  process  (A^"  +  e  — >  A  +  A*)  is  also  included,  that  is  the 
fastest  mechanism  of  bulk  recombination  for  molecules,  with  a  rate  of  kr 3  =  1.0  x  10_13TfT°-5m3/s.  It  can 
be  seen  that  the  plasmoid  translation  and  entrainment  is  expected  to  be  dominated  by  recombination  at 
temperatures  below  ~5  eV  and  by  the  electron  impact  ionization  at  higher  temperatures.  The  latter  one  is 
related  to  a  relatively  low  ionization  threshold  of  15.6  eV  as  compared  to  the  average  energy  in  translational 
and  internal  modes  of  molecular  nitrogen.  The  SCX  rate  for  nitrogen  is  significantly  lower  than  that  of  the 
other  two  processes  over  the  whole  temperature  range  of  interest,  so  the  use  of  nitrogen  in  FRC  thruster 
may  be  problematic. 


Figure  2.  Plasma  reaction  rates  in  helium  (left)  and  molecular  nitrogen  (right). 

High  thrust  levels  can  be  achieved  using  heavy  particles,  so  xenon  is  often  used  in  thrusters.  However, 
because  it’s  ionization  energy  is  low,  the  Eli  rate  for  xenon  is  about  two  orders  of  magnitude  larger  than 
that  of  SCX  over  the  whole  temperature  range  (see  Fig.  2,  left).  Therefore,  the  efficiency  of  xenon-based 
FRC  thruster  is  expected  to  be  low. 

Species  with  high  ionization  energy  have  low  Eli  rate.  One  of  such  species  is  neon,  whose  ionization 
energy  is  equal  to  21.6  eV.  The  charge  exchange  and  ionization  rates  for  neon  are  plotted  in  Fig.  3  (right). 
The  recombination  rate  is  not  shown  as  it  was  found,  similar  to  helium,  to  be  much  lower  than  the  other 
two  rates.  The  results  show  that  the  charge  exchange  reaction  is  expected  to  occur  more  often  than  the 
ionization,  and  for  plasmoid  temperatures  on  the  order  of  10  eV  or  lower  the  entrainment  efficiency  may  be 
sufficiently  high.  Therefore,  neon  will  be  used  in  all  subsequent  calculations. 
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Figure  3.  Plasma  reaction  rates  in  xenon  (left)  and  neon  (right). 

III.  Adiabatic  relaxation  of  neon  plasma 

When  cold  stagnant  neutral  atoms  collide  with  the  fast  moving  plasma  ions  and  electrons,  the  relaxation 
process  that  ensues  depends  on  many  factors  such  as  plasma  temperature,  concentration  of  reacting  species 
and  the  relationship  between  different  reaction  rates.  Let  us  consider  the  effect  of  these  factors  on  the  time 
evolution  of  plasma  and  neutral  properties  in  an  adiabatic  heat  bath.  Although  the  spatially  homogeneous 
heat  bath  can  not  accurately  predict  the  multi-dimensional  evolution  of  plasma  affected  by  the  neutral 
flow,  it  nevertheless  captures  the  main  features  of  the  relaxation  process,  and  thus  helps  to  understand  the 
relative  importance  of  various  factors  at  play.  To  better  model  the  neutral  entrainment  process,  it  is  assumed 
that  initially  the  neutral  species  and  the  charged  species  have  relative  velocity  U .  The  dependence  of  the 
relaxation  process  characteristics  on  the  magnitude  of  U  is  illustrated  in  Fig.  4.  Unless  otherwise  noted,  the 
initial  plasma  temperature  is  set  to  10  eV,  the  neutral  temperature  is  300  K,  and  both  neutral  and  plasma 
density  are  1018  m~3.  The  calculations  are  performed  for  the  two  different  values  of  U,  20  krn/s  and  30 
km/s. 

Consider  first  the  case  of  20  km/s.  In  the  absence  of  relative  collision  velocity,  the  ions  and  neutrals  are 
expected  to  quickly  relax  to  0.5  *  ( Tneut  +Tion )  «  5.01  eV.  However,  as  the  kinetic  energy  of  relative  motion 
of  charged  and  neutral  particles  is  converted  into  thermal  energy,  the  temperature  of  all  species  increases, 
and  heavy  particles  achieve  thermal  equilibrium  at  the  temperature  of  about  12  eV  (see  Fig.  4,  left).  The 
thermal  relaxation  of  electrons  due  to  elastic  collisions  with  neutrals  and  Coulomb  collisions  with  ions  is 
fairly  slow,  and  the  change  in  electron  temperature  is  primarily  related  to  the  electron  impact  ionization 
reactions,  which  ,  being  strongly  endothermic,  significantly  decreases  electron  temperature. 

For  a  higher  relative  velocity  between  charged  and  neutral  species  (30  km/s),  due  to  higher  values  of  the 
kinetic  energy  of  relative  motion,  the  final  temperature  of  heavy  species  is  over  20  eV.  The  electron  temper¬ 
ature  changes  very  little  compared  to  the  20  km/s  case.  This  is  again  related  to  the  slow  thermal  relaxation 
of  electrons,  and  small  change  in  ionization  rate  as  the  relative  electron/neutral  velocity  is  small  compared 
to  the  thermal  velocities  of  electrons,  and  the  latter  is  the  main  contributor  to  ionization.  The  impact  of 
relative  velocity  on  electron  temperature  becomes  visible  only  after  100  /is.  The  important  conclusion  that 
can  be  drawn  from  these  results  is  that  for  any  flow  velocity  (and  also  for  any  plasma  temperature)  there  is  a 
strong  thermal  non-equilibrium  between  neutrals,  electrons,  and  ions.  Such  a  non-equilibrium  is  expected  to 
occur  in  the  actual  neutral  entrainment  of  an  FRC  thruster,  and  therefore  a  kinetic  approach  is  best  suited 
to  model  such  an  entrainment. 

The  number  density  of  ions  very  weakly  depends  on  the  flow  speed,  (see  Fig.  4, right)  because  number 
of  ions  changes  only  due  to  Eli  reactions,  and,  as  stated  earlier,  their  rates  are  not  much  affected  by  the 
flow  speed.  The  initial  electron  temperature,  however,  is  extremely  important,  as  shown  in  Fig.  4  (right) 
for  U  =  20  km/s.  The  electron  temperature  impacts  the  ionization  rate,  and  when  the  initial  plasma 
temperature  is  5  eV,  less  than  3%  of  the  neutral  gas  is  ionized  in  about  30  fis.  To  study  the  dependence  of 
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the  relaxation  process  on  the  initial  density  of  the  neutrals,  the  calculation  was  performed  with  the  initial 
number  density  of  the  neutrals  equal  to  3xl018  m-3.  The  number  of  charge  exchange  reactions  for  a  plasma 
temperature  of  10  eV  and  U  =  20  km/s  was  found  to  linearly  increase  with  neutral  density.  The  dependence 
of  the  number  of  ionization  reactions  on  neutral  density  is  weaker  than  linear,  as  can  be  seen  in  Fig.  4 
(right).  When  neutral  density  increases  by  a  factor  three,  the  ion  density  increase,  primarily  affected  by  the 
ionization  process,  is  roughly  a  factor  of  three  only  in  the  first  few  microseconds.  At  1000  /is,  the  ion  density 
change  is  only  25%  higher  than  the  baseline  case.  This  is  primarily  related  to  the  depletion  of  high  energy 
electrons,  as  will  be  discussed  below. 


Time,  |rs  Time,  its 

Figure  4.  Impact  of  plasma  parameters  on  the  temperature  of  different  species  (left)  and  ion  number  density 
(right) 

Generally,  the  electron  impact  ionization  reactions  results  in  energy  transfer  from  the  translational  motion 
of  reacting  electrons  and  neutrals  to  the  potential  energy  of  newly  formed  electron  and  ion.  The  reaction 
therefore  decreases  plasma  temperature,  and  in  particular,  depletes  the  high  velocity  tail  of  electrons.  This, 
in  turn,  will  slow  down  the  electron  impact  ionization  rate  (but  have  only  limited  impact  on  charge  exchange 
reactions).  The  main  process  that  leads  to  the  equilibration  of  the  electron  velocity  distribution  function 
is  Coulomb  collisions  between  electrons.  When  the  Coulomb  collisions  are  turned  off,  the  high  energy  tail 
is  quickly  depleted  by  ionization,  as  illustrated  in  Fig.  5  (left).  There  is  also  noticeable  impact  of  Coulomb 
collisions  on  the  evolution  of  plasma  density,  shown  in  Fig.  5  (right)  for  U  =  20  km/s  and  Tpini  =  10  eV. 
Note  that  such  a  non-equilibrium  in  velocity  distribution  is  another  indication  that  a  kinetic  approach  is 
necessary  to  accurately  model  neutral  entrainment  in  an  FRC  thruster.  For  a  multi-dimensional  neutral 
entrainment,  when  a  translated  plasmoid  passes  through  neutral  gas,  the  effects  mentioned  in  this  section, 
such  as  non-equilibrium,  fast  relaxation  rates  for  high  neutral  density,  strong  influence  of  plasma  temperature, 
the  increase  of  ionization  rate  due  to  Coulomb  collisions,  etc,  are  expected  to  be  even  more  pronounced,  as 
the  plasma  particles  collide  with  new,  undisturbed  neutral  gas  due  to  the  plasmoid  translation. 

IV.  Numerical  method  and  flow  conditions  for  FRC  entrainment 

The  computational  tool  that  was  used  in  this  work  to  model  the  process  of  neutral  entrainment  of  an  FRC 
plasmoid  is  Celeste3D.  Celeste3D  is  an  implicit  three-dimensional  Particle-In-Cell  (PIC)  code  that  solves  the 
full  set  of  Maxwell- Vlasov  equation  and  has  been  extensively  used  in  the  past  for  various  astrophysical  and 
laboratory  plasma  problems.12, 13  The  implicit  moment  formulation  of  the  PIC  method  implemented  in 
Celeste  results  in  highly  efficient  simulations  based  on  ion  length  and  time  scales,  and  not  electron  scales 
as  explicit  methods  do,  while  retaining  the  kinetic  effects  of  both  the  electrons  and  the  ions.  An  explicit 
simulation  requires  the  time  step  to  be  At  <  2/copei  and  the  spatial  cell  size  to  be  Ax  <  /Xe  in  order  to 
avoid  the  finite  grid  instability.  Here,  uipe  is  the  electron  plasma  frequency,  and  Ae  is  the  electron  Debye 
length.  In  an  implicit  simulation,  these  requirements  are  replaced  by  an  accuracy  condition  related  to  the 
conservation  of  energy,  At<Deltax  where  ce  is  the  electron  thermal  speed. 

Original  Celeste3D  models  the  evolution  of  ion  and  electron  populations  as  described  by  the  Vlasov 
equation,  coupled  with  the  solution  of  full  Maxwell  equations;  it  does  not  simulate  neutral  particles.  To 
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Figure  5.  Depletion  of  electron  velocity  distribution  function  due  to  ionization  when  Coulomb  collisions  are 
turned  off  (left)  and  the  effect  of  Coulomb  collisions  on  plasma  density  evolution  (right). 


study  the  neutral  entrainment  process,  Celeste3D  was  extended  to  include  neutral  transport  and  collisional 
relaxation.  A  Direct  Simulation  Monte  Carlo  based  capability  has  been  added  to  Celeste  that  includes 
neutral  transport  and  collisions.  The  following  collision  processes  are  included  in  this  work:  neutral-neutral 
collisions  (modeled  according  to  the  VHS  model15);  the  charge  exchange  reactions  (the  cross  section  is  given 
earlier  in  this  paper);  neutral-ion  elastic  collisions  (according  to  Ref.,10  the  cross  section  is  twice  the  charge 
exchange  cross  section);  the  electron  impact  ionization  (the  cross  section  is  also  given  earlier).  The  hard 
sphere  after-collision  scattering  is  assumed  for  all  these  processes,  with  the  exception  of  charge  exchange 
reactions,  for  which  the  velocities  of  neutrals  and  ions  are  swapped.  As  all  species  have  different  weights,  a 
weighting  scheme16  is  applied.  The  majorant  collision  frequency17  of  the  DSMC  method  is  used  to  model 
the  collision  process  in  cells.  In  addition  to  the  neutral  capability,  a  Coulomb  collision  module  has  been 
added  to  Celeste,  based  on  a  particle- weights  scheme  of  Ref.18  (see  previous  section  for  more  detail). 

Although  the  actual  entrainment  geometry  is  three-dimensional,  as  a  nearly  axisymmetric  plasmoid 
interacts  with  neutral  gas  supplied  from  two  or  more  azimuthal  injectors,  three-dimensional  PIC  modeling 
is  extremely  time  consuming,  and  in  this  work  it  was  replaced  by  a  2D  (planar  flow).  The  result  is  a  much 
larger  number  of  particles  per  cell  and  thus  smaller  statistical  scatter.  The  grid  was  40x80,  and  the  number 
of  simulated  ions,  electrons,  and  neutrals  per  cell  was  approximately  64,  100,  and  512,  respectively.  Since 
the  cost  of  an  implicit  simulation  is  a  direct  function  of  the  electron  mass  (see  above,  the  computational 
time  is  inversely  proportional  to  the  square  root  of  the  electron  mass  me ),  the  simulation  efficiency  may  be 
further  enhanced  through  the  introduction  of  the  weighted  electron  mass,  m'e  =  Wme,  where  the  constant 
W  >  1.  The  ion-to-electron  mass  ratio  has  to  be  high  enough  to  preserve  the  kinetic  effects,  and  the  value  of 
rrii/m'e  on  the  order  of  100  is  usually  sufficient.12  In  this  work,  mi/m'e  =  100  is  used.  The  initial  conditions 
used  to  set  an  equilibrium  plasmoid  is  the  a  Schmid-Burgk  equilibrium14  with  the  average  plasma  density  of 
1018  nr-3  and  two  plasma  temperature  of  5  and  10  eV. 

To  simplify  the  inflow  boundary  conditions,  the  neutral-plasma  interaction  is  examined  in  the  reference 
frame  of  the  plasmoid,  so  that  the  neutral  gas  is  injected  into  the  computational  domain  from  the  left 
boundary,  and  then  passes  through  the  domain  from  left  to  right.  The  properties  of  the  injected  neutral  gas 
are  changed  to  study  their  impact  on  the  neutral-plasma  interaction.  In  these  computations,  both  periodic 
(10  eV  plasmoid)  and  open  (5  eV  plasmoid)  boundary  conditions  for  plasma  and  open  boundary  conditions 
for  neutrals  were  imposed  at  the  left  and  right  boundaries,  and  a  conducting  wall  with  specular  reflection 
was  set  at  the  top  and  bottom  boundaries.  Note  that  computations  were  also  performed  with  no  neutral 
inflow  to  analyze  the  temporal  evolution  of  the  plasmoid;  it  was  observed  that  the  equilibrium  is  maintained 
for  at  least  10,000  ion  plasma  periods.  This  amounts  to  over  30  fxs  for  the  conditions  under  consideration, 
and  is  longer  than  the  plasmoid-neutral  interaction  time. 
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V.  Interaction  of  FRC  plasmoid  with  neutral  gas 


The  first  set  of  computations  is  conducted  for  the  average  plasma  density  inside  the  plasmoid  of  1018  molecule/: 
and  plasma  temperature  of  5  eV,  with  the  neutral  density  of  1018  molecule/nr3  and  a  relative  velocity  Urei 
of  20  krn/s  between  the  plasmoid  and  neutrals.  As  mentioned  earlier,  the  plasmoid-neutrals  interaction  is 
modeled  in  the  reference  frame  of  the  plasmoid,  so  that  the  plasmoid  is  initially  not  moving,  and  neutral  gas 
is  introduced  into  the  computational  domain  with  the  bulk  velocity  of  Urei.  Time  evolution  of  neutral  gas 
density  is  shown  in  Fig.  6  (left).  The  gas  is  injected  from  the  left  boundary,  and  the  density  is  normalized  by 
the  free  stream  value.  At  the  time  moment  of  8  /is  after  the  beginning  of  the  plasmoid/neutrals  interaction, 
the  neutrals  advanced  past  the  center  of  the  computational  domain  into  the  region  where  the  plasma  density 
is  maximum.  Although  there  is  some  small  decrease  in  neutral  density  due  to  electron  impact  ionization 
reactions  as  well  as  some  backward  scattering,  this  effect  is  fairly  small  at  this  point.  However,  at  12  fj, s  it 
becomes  more  obvious,  with  the  loss  in  neutral  density  reaching  about  20%  in  the  central  part  of  the  flow. 

As  the  neutrals  reach  the  outflow  boundary  at  time  16  /is,  this  loss  amounts  to  about  30%.  Again,  the  loss 
is  both  due  to  elastic  and  inelastic  interactions  between  neutral  atoms  and  charged  particles. 

The  field  of  ion  density  normalized  by  1018  molecule/m3  is  shown  in  Fig.  6  (right)  at  the  same  time 
moments.  Note  that  the  electron  density  is  nearly  identical,  except  for  some  statistical  scatter  inherent  in 
the  approach.  After  the  first  8  /xs,  the  impact  of  neutrals  on  the  plasmoid  is  negligible  as  the  ion  density 
is  almost  symmetric  with  respect  to  X  =  0.15  m.  At  12  fi s,  the  impact  of  neutrals  on  the  plasmoid  is 
manifested  through  a  moderate,  on  the  order  of  5%,  increase  in  plasma  density  in  the  region  of  maximum 
plasma  density  as  a  result  of  electron  impact  ionization,  and  a  significant  elongation  of  the  plasmoid  due  to 
the  momentum  transfer  between  plasma  and  neutrals  in  charge  exchange  and  elastic  collisions.  This  effect 
is  even  more  pronounced  at  16  /is,  when  the  plasma  density  in  the  center  decreases  as  the  result  of  mass 
transfer  in  the  direction  of  the  neutral  flow,  Open  boundary  conditions  at  the  left  and  right  boundaries 
are  used  for  these  computations,  and  a  significant  loss  of  plasma  particles  though  the  outflow  boundary  is 
observed  at  16  /its, 


Consider  now  the  results  of  the  interaction  between  a  10  eV  plasmoid  and  neutral  gas.  Figure  7  (left) 
presents  an  instantaneous  snapshot  of  ion  density  at  a  time  of  15  p,s  after  the  neutral  injection  started 
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for  all  cases  except  U  =  30  km/s,  for  which  it  is  10  /is  after  the  start  of  injection  (this  is  the  time  that 
the  moving  neutrals  need  to  cross  the  computational  domain).  The  ion  number  density  presented  in  that 
figure  is  normalized  by  1018.  For  comparison,  the  result  for  a  no-neutrals,  plasma-only  case  is  also  shown 
(top  quarter).  Note  that  the  no-neutral  field  practically  does  not  change  with  time,  with  the  exception  of 
small  changes  related  to  the  statistical  scatter.  The  baseline  case  shown  in  the  second  quarter  represents 
the  free  stream  neutral  density  of  1018  and  velocity  of  20  km/s,  and  the  third  and  forth  cases  consider  an 
increased  neutral  velocity  and  density,  respectively.  For  the  baseline  case,  the  neutral-plasma  interaction 
results  in  modest,  slightly  more  than  5%,  increase  in  the  number  of  ions  in  the  region  where  the  ion  density 
is  near  its  maximum.  This  is  obviously  the  result  of  electron  impact  ionization  reactions.  It  is  important 
to  note  that  there  is  a  clear  translation  of  the  plasmoid  from  the  center  of  the  computational  domain  to 
the  right.  The  shift,  which  amounts  to  almost  2  cm,  is  due  to  the  momentum  transfer  from  the  moving 
neutrals  to  the  plasma  particles.  The  average  velocity  of  initially  stationary  plasmoid  is  about  3  km/s.  Since 
the  computations  are  conducted  in  the  frame  reference  of  the  plasmoid,  in  application  to  the  actual  FRC 
entrainment  configuration  that  means  that  the  momentum  is  in  fact  transferred  from  the  moving  plasmoid 
to  the  neutral  gas.  For  a  higher  relative  velocity  between  the  plasmoid  and  neutrals,  U  =  30  km/s,  the 
plasma-neutral  interaction  is  clearly  weaker.  The  plasmoid  is  translated  only  at  about  1  cm,  and  the  number 
of  neon  ions  is  increased  only  by  about  3%.  The  reason  for  this  is  a  shorter  interaction  time  due  to  a  higher 
relative  velocity  between  neutral  and  plasma  particles.  When  the  free  stream  neutral  density  is  increased 
by  a  factor  of  three,  the  momentum  transfer  essentially  triples,  as  does  the  number  of  ionization  reactions. 
Note  also  that  the  shape  of  the  plasmoid  becomes  more  elongated  for  stronger  neutral-plasma  interaction. 
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Figure  7.  Ion  density  (left)  and  temperature  (right)  for  different  free  stream  conditions. 


The  ion  temperature  increases  for  the  baseline  case  by  about  1.6  eV  in  most  of  the  computational  domain 
as  compared  to  the  plasma-only  case,  as  illustrated  in  Fig.  7  right).  Here  and  below,  the  temperature  is  in  eV. 
The  increase  is  much  more  pronounced,  over  10  eV,  for  the  U  =  30  km/s  case.  This  is  related  to  the  transfer 
of  energy  from  the  bulk  neutral  flow  motion  to  the  thermal  motion  of  ions  and  neutrals,  and  is  consistent 
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with  the  conclusions  of  the  Heat  Bath  section.  Note  that  for  U  =  30  km/s  the  maximum  ion  temperature 
is  observed  to  the  right  of  the  maximum  ion  density,  where  more  ions  have  collided  with  electrons.  For 
the  elevated  neutral  density,  there  is  a  minimum  ion  temperature  in  this  region,  related  primarily  to  the 
large  number  of  endothermic  ionization  reactions  (for  the  baseline  case,  there  is  a  temperature  minimum 
in  this  region,  too,  although  less  pronounced).  The  number  density  of  neutral  particles,  normalized  by  its 
free  stream  value,  is  plotted  in  Fig.  8  (left).  It  can  be  seen  that  the  neutral  flow,  moving  from  left  to 
right,  almost  reached  the  right  boundary  of  the  computational  domain.  Most  importantly,  the  interaction 
of  neutral  particles  with  ions  and  electrons  depleted  their  population  in  the  central  region  by  almost  50% 
for  the  baseline  and  elevated  density  cases,  and  by  over  30%  for  U  =  30  km/s.  Recalling  that  the  absolute 
change  in  ion  density  is  much  smaller  than  the  change  in  neutral  density,  it  becomes  clear  that  the  latter  is 
related  to  charge  exchange  reactions.  The  neutral  particles  loose  there  momentum  in  X  direction  after  charge 
exchange  reactions,  and  thus  do  not  reach  the  right  boundary.  The  temperature  field  of  neutral  particles, 
shown  in  Fig.  8  (right),  indicates  the  region  where  neutrals  mostly  travel  after  they  are  created  as  a  result  of 
charge  exchange  reactions.  It  is  the  region  of  elevated  temperature  (relative  velocity  between  these  neutrals 
and  free  stream  neutrals  is  high)  observed  mostly  to  the  left  of  the  center  of  the  plasmoid,  and  away  from  the 
top  and  bottom  conducting  walls.  Therefore,  in  the  actual  entrainment  the  neutrals  are  expected  to  travel 
with  the  plasmoid.  Obviously,  this  is  exactly  the  effect  that  the  entrainment  is  used  for. 


Figure  8.  Neutral  density  (left)  and  temperature  (right)  for  different  free  stream  conditions. 


VI.  Conclusions 

This  work  is  a  step  toward  accurate  modeling  of  a  high-power  electric  propulsion  thruster  based  on 
a  high-speed  neutral  entrainment  of  a  field  reversed  configuration  plasmoid.  A  number  of  chemical  and 
physical  processes  are  important  in  the  interaction  between  the  plasmoid  and  neutral  gas,  and  two  of  them 
are  considered  here,  charge  exchange  and  electron  impact  ionization  reactions.  Both  processes  are  expected 
to  significantly  affect  the  performance  of  an  FRC  thruster,  with  the  former  increasing  thrust  force  and 
improving  total  thruster  efficiency  and  the  latter  is  likely  decreasing  them.  Comparison  of  ionization  and 
charge  exchange  reaction  rates  presented  in  this  work  indicates  that  the  use  of  nitrogen  and  xenon  may  be 
problematic,  while  neon  appears  to  be  a  fairly  good  propellant. 

The  computation  effort  presented  in  this  work  covered  thermal  relaxation  in  an  adiabatic  heat  bath  with 
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neutral  and  plasma  properties  close  to  those  of  a  typical  FRC  thruster,  and  a  two-dimensional  interaction  of 
a  1018  molecule/m3,  5  to  10  eV  FRC  plasmoid,  initially  at  a  Schmid-Burgk  equilibrium,  and  neutral  gas  with 
comparable  density.  The  0D  computations  showed  that  the  relaxation  process  proceeds  under  conditions  of 
strong  thermal  and  chemical  non-equilibrium;  ion,  electron,  and  neutral  temperatures  strongly  differ,  and 
the  electron  distribution  function  is  non-Maxwellian.  This  indicates  that  a  kinetic  approach  has  to  be  used 
to  model  neutral  entrainment  in  FRC  thrusters.  Strong  impact  of  electron  temperature  on  plasma  density  is 
shown  in  heat  bath  relaxation,  mostly  related  to  ionization  reactions.  Increase  in  relative  velocity  between 
plasma  and  neutrals  strongly  affects  ion  and  neutral  temperature.  Electron  impact  ionization  was  found  to 
deplete  high  velocity  tail  of  the  electron  distribution  function,  and  modeling  of  Coulomb  collisions  between 
electrons  is  desirable  to  properly  account  for  that  depletion. 

Two-dimensional  modeling  of  a  an  FRC  plasmoid  /  neutral  gas  interaction  is  conducted  with  an  implicit 
PIC  code  Celeste3D,  extended  in  this  work  to  include  neutral  transport,  plasma-neutral  and  neutral-neutral 
collisions  and  Coulomb  collisions.  Neon  propellant  was  used  with  a  baseline  plasma  and  neutral  densities 
of  1018  nr-3,  plasma  temperatures  of  5  and  10  eV,  and  relative  velocities  of  20  and  30  km/s.  For  both 
plasma  temperatures,  the  results  show  strong  entrainment  of  neutral  particles  by  a  translated  plasmoid  as 
a  result  of  charge  exchange  reactions  between  slow  neutrals  and  fast  moving  ions,  and  a  modest  increase  in 
plasma  density  due  to  electron  impact  ionization.  The  increase  in  neutral  density  appears  highly  beneficial 
in  terms  thruster  efficiency  as  it  proportionally  increases  the  number  of  entrained  neutrals.  The  increase 
in  plasmoid  velocity  decreases  that  number  due  to  a  shorter  plasma-neutral  interaction  time.  Not  that  the 
present  computations  have  not  included  the  electronic  excitation  and  related  processes,  the  effect  of  those  is 
estimated  to  be  similar  to  the  electron  impact  ionization  for  plasmoid  temperatures  between  5  and  10  eV, 
and  will  be  included  in  the  future. 
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